Cite this Article
Su Guifeng, Li Xiaowen, Zhang Xiaobing, Zhang Yi. Finite density scaling laws of condensation phase transition in zero-range processes on scale-free networks. Chinese Physics B, 2020, 29(8): 088904  
Permissions
Finite density scaling laws of condensation phase transition in zero-range processes on scale-free networks
Su Guifeng1, Li Xiaowen1, Zhang Xiaobing2, Zhang Yi1, †
Department of Physics, Shanghai Normal University, Shanghai 200234, China
School of Physics, Nankai University, Tianjin 300071, China

 

† Corresponding author. E-mail: yizhang@shnu.edu.cn

Project supported by the National Natural Science Foundation of China (Grant No. 11505115).

Abstract

The dynamics of zero-range processes on complex networks is expected to be influenced by the topological structure of underlying networks. A real space complete condensation phase transition in the stationary state may occur. We study the finite density effects of the condensation transition in both the stationary and dynamical zero-range processes on scale-free networks. By means of grand canonical ensemble method, we predict analytically the scaling laws of the average occupation number with respect to the finite density for the steady state. We further explore the relaxation dynamics of the condensation phase transition. By applying the hierarchical evolution and scaling ansatz, a scaling law for the relaxation dynamics is predicted. Monte Carlo simulations are performed and the predicted density scaling laws are nicely validated.

PACS: ;89.75.Hc;;05.20.-y;;02.50.Ey;
Keyword:;finite density scaling;;zero-range processes (ZRP);;scale-free networks;
1. Introduction

Condensation phase transitions are abundant phenomena in nature and observed in various physical circumstances and stochastic processes, for instance, the Bose–Einstein condensation (BEC) in cold dilute atomic gases,[1] jamming in traffic flow[2] and granular flow,[3,4] condensation of links in complex networks.[5–7] Although the equilibrium condensation phase transitions can be well described in the framework of equilibrium statistical mechanics,[8] building a consistent theory of nonequilibrium (condensation) phase transitions remains a big challenge. How such transitions arise and what dynamical nature they may have are particularly intriguing questions to be answered. The long lasting efforts have been taken to investigation of interacting particle systems.[9] One of the important categories is so-called driven diffusive systems.[10] These systems are driven away from equilibrium by external forces, and may evolve toward a nonequilibrium steady state. In addition, they exhibit non-trivial behavior such as phase separation[11] or phase transitions (even in one spatial dimension).[12]

Among the above-mentioned systems, the zero-range processes (ZRPs), introduced by Spitzer[13] roughly half century ago, has attracted considerable attention due to the fact that it is one of a few rare examples of analytically tractable nonequilibrium models. For recent reviews, see e.g., Refs. [14,15] and the references therein. The prototype of ZRP is a stochastic particle system defined on one-dimensional regular lattice of given number of sites and particles hop onto the nearest neighbors (NNs) with some transition rate, say, u(n). This rate function depends only on the occupation number n of its current position, rather than the other sites. This defines the zero-range property of the process. If u(n) is decreasing in n, then an effective attraction between particles is presented. As a result, a condensation transition will occur if ρ exceeds some critical density ρc, and a macroscopically finite fraction of particles condense on a single site. The condensation phase transitions in ZRPs have been studied based on periodic lattices in Refs. [16,19].

From the perspective of complex networks, the underlying periodic lattice of ZRPs is simply regular fully connected network. However, many natural and manmade networks are self-organized as the scale-free (SF) networks (see e.g., Refs. [20–22] and the references therein). The SF networks are strongly inhomogeneous and highly clustering in architecture.[23] They also typically possess a power-law degree distribution, , where k is the degree (i.e. number of links) of a node and the constant γ is degree exponent. It has been demonstrated that these significant features of SF networks are important factors affecting a multitude of critical phenomena,[24] and as well as dynamical processes[25–29] defined on them. It was found that when the power exponent δ in the hopping rate u(n) = nδ is smaller than a critical value δc, that is, when δ < δc = 1/(γ – 1), a finite fraction of particles are condensed into the node with the maximum degree, or the hub.[27] Such condensation phase transitions are driven by the interplay of the (effective) interaction between particles and the topological structure disorder of the underlying SF networks. However, to the author’ knowledge, so far the analysis of the possible effects of the finite density are still absent. In fact, one naturally expects the density plays an important role in condensation transition of ZRPs on SF networks.

In this paper we then focus on the finite density aspect of the condensation transition of the static and dynamical ZRP defined on SF networks. As we will show in the following, based on the grand canonical ensemble (GCE) approach, we predict a novel scaling law in the condensation phase transition with respect to (w.r.t.) the finite density. Another important issue we addressed in this study is the relaxation dynamics in the condensed phase transition. A hierarchical evolution occurs in relaxation dynamics, i.e., the relaxation proceeds from small degree nodes to larger and larger degree nodes. We predict that there exist the scaling laws of the average occupation number and the inverse participation ratio (IPR) w.r.t. the finite density in both steady state and relaxation dynamics. To verify our analytical predictions, Monte Carlo simulations are fulfilled for steady state and the relaxation dynamics, respectively. These predicted finite density scaling laws are nicely confirmed.

The rest of the paper is organized as follows: in Section 2, we present the ZRP model on SF network, and the corresponding GCE theoretical formalism. In Section 3, we first argue the existence of the condensation phase transition, and then predicted the finite density scaling laws for such a transition. Monte Carlo simulations are followed and fully confirmed previous analytical results. In Section 4, we explore the relaxation dynamics both analytically and numerically. Predictions on the finite density scaling law for relaxation dynamics are presented and simulated. Finally in Section 5, we conclude the paper with the summary of our findings and some remarks.

2. ZRP on SF network
2.1. The model

Let us now consider the ZRP with N = ρ L particles hopping on an SF network, with ρ being the density of the particles, and L the number of nodes of the underlying SF networks. One has to keep in mind that in current situation, the density ρ is not fixed but variable. The degree (probability) distribution of the SF networks is characterized by power-law distribution

where c is some appropriate normalization constant. For any finite network, the minimum degree kmin sets up the smallest possible degree and is a constant of order one, . Recall that the node with the maximum degree is called hub, and its degree kmax satisfies[30]

such that kmax ≈ kmin L1/(γ – 1). In general, one has 1 ≤ kmin ≪ kmax. We use wi to represent the weight of the i-th node of the network and wi = ki/Σi ki ≡ ki/Σk, where a shorthand notation Σk is used for convenience (the same below). The average degree of SF networks is defined as

To keep finite, the degree exponent γ > 2 is required.

One of the important properties of the ZRP stationary state probability , i.e., the probability of finding particles in a microscopic configuration {ni}, is that takes the factorized form,[31,32] i.e.,

in which the microscopic configuration {ni} = {n1, n2, …, nL}, ni is the number of particles on the i-th node of the network. The normalization factor ZL,N plays the role of partition function, and can be computed by summing the product over all microscopic configurations {ni},

The δ function introduced in Eq. (5) guarantees the conservation of total number of particles by the dynamics, i.e., Σn = N. The factors appeared in the product form f(nj)’s are the single-site (the j-th node) statistical weights, and can be expressed in terms of hopping rate functions u(nj),

where Wj is the single particle stationary state probability distribution, associated with the jumping probability Ti → j. It can be written as

The subscript i → j in jumping probability means that particles hop from the i-th (departure) node to the j-th (arrival) node. The detailed balance equations fully determine the Wj’s:

i.e., the outgoing probability flux exactly matches the incoming probability flux. The solution to Eq. (8) is Wi = ki/Σk = wi,[31] i.e., the stationary state probability distribution just equals the weight of the node: the more weights a node has, the more frequent it will be visited by a random hopper. For the ZRP on the SF networks, the hopping rate function of particles out of the j-th node takes the form

where nj represents the number of particles located in the current (j-th) node, and δ is a global parameter controlling the interaction strength between particles. For particles with (effective) attractive interaction, u(n) grows sub-linearly w.r.t the number of particles n, while for (effective) repulsive interaction, u(n) grows super-linearly w.r.t n. As a result, the f(nj) turns out to be

2.2. Grand canonical ensemble formalism

We now present the general GCE approach that deals with the condensation phase transition in ZRPs.[16] The grand canonical partition function is defined by

where z is the fugacity. Physically it is the average hopping rate. Applying Eq. (5) to Eq. (11) one obtains

where is

The average occupation number of the i-th node mi reads

The self-consistent equation ρ = Σ m/L determines the fugacity. By introducing the effective fugacity xi,

and with help of Eqs. (4), (10), and (13), mi can be written as

in which

We will apply the above GEC approach to the ZRP on SF networks, and this is the major task in the following sections.

3. Finite density scaling law of the steady state
3.1. Analytical results

The so-called complete condensation phase transition in ZRPs on SF network with a given density ρ has been revealed in a recent study.[27] The structural inhomogeneity of SF networks plays an important role. To summarize, when the complete condensation occurs, almost a whole fraction of particles are condensed onto a few high-degree nodes. As we will see in the following, a variable finite density has some nontrivial effects on the condensation phase transitions.

Although the inter-particle interactions are determined by the parameter δ in hopping rate function, only the case of 0 < δ < 1 is what we are concerned in this paper. When δ = 0, the model is equivalent to the disordered ZRP with n-independent hopping rate functions, and has already been studied in Refs. [32–34]. The analytical solution in this case is simply mi = xi/(1 – xi). When δ = 1, the hopping rate function becomes u(n) = n, which is proportional to the number of particles. Physically this corresponds to fully independent motion of particles, and the system consists of N non-interacting random particles. One may solve the occupation number distribution mi(z) = zki,[27] as is expected.

We then pay our attention to 0 < δ < 1. To solve the average occupation number mi in this general case, one has to resort to the approximate expression of in Eq. (17), since in this case the analytically closed form does not exist. The saddle-point approximations can be borrowed to finish the task for the series in Eq. (17). For large x, it can be shown[27] that

For small x, the series in Eq. (17) can be approximated by summing over a few lowest order terms. As a result, one obtains mi ≃ xi for xi ≪ 1, and for xi ≳ 1.

The hub owns the most links to another nodes, its fugacity is . The effective fugacities of the rest of nodes are , where kc is the crossover degree, and defined as

The minimum effective fugacity is obviously xmin = kmin/kc. The kc is determined by the self-consistency condition ρ = Σ m/L and the possible form of the minimum effective fugacity xmin. Let us assume firstly xmin ≳ 1, then the self-consistency condition reads

The in the above equation represents the average of (…), is just

by definition. In thermodynamic limit L → ∞, the summation (over all nodes) in Eq. (21) can be replaced by integral, hence

By means of Eq. (20), the crossover degree kc w.r.t. the density can be obtained, i.e.,

Note that kc is density dependent and scales as kc ∼ ρ−δ. The average occupation number mk on a node with degree k > kc for steady state becomes

On the hub, it scales as

where the last relationship comes from the fact that kmax ≈ kmin L1/(γ – 1) ∼ Lδc for the SF network. Since δc/δ < 1, the sub-linear increasing of mhub on L suggests that no condensation transition occurs, the whole system is in a single non-condensation phase. Physically, this is due to the existence of the lower bound degree kmin. When kc < kmin, the condensation transition is fully suppressed and only a single {fluid-like} phase exists for δ > δc.

On the other hand, if we assume xmin ≪ 1, or equivalently kc ≫ kmin, the density can be decomposed into a fluid part ρf and a condensed part ρc:

with

In the thermodynamic limit L → ∞, the summands in ρc and ρf read, respectively,

The integral is finite such that a condensation phase transition may occur if ρf vanishes as . As a result,

The self-consistency condition requires the integral in Eq. (30) divergent as L → ∞. This is to require the exponent part in the above integral, 1/δ – γ > – 1, such that δ ≤ δc. One sees that kc is the degree at which the fluid-like phase crossovers to the condensation phase, and scales as kc ∼ (ln kmax/ρ)δc for δ = δc. For δ < δc, it scales as

For an SF network with finite size L, Eq. (31) is reduced to kc ∼ Lδc – δ/ρδ for δ < δc. Note that the presence of the factor ρ−δ implies a decreasing of kc as the density increases. Further, the average occupation number scales as mk ∼ ρδck/(ln kmax)δc (k < kc), and mk ∼ ρ k1/δc/ln kmax (k > kc) for δ = δc.

For δ < δc, we have

For a finite number of nodes L, Eq. (32) becomes

for the fluid-like phase (k < kc) and the condensation phase (k > kc), respectively. In contrast to the case with a fixed density, now the density ρ comes in and has effects on the average occupation number and hence the condensation phase transition.

We summarize our predicted results in Table 1 for readability. The finite density scaling relations are explicitly presented.

Table 1.

The scaling relations of the crossover degree kc and the average degree occupation number mk w.r.t. the finite density ρ and the parameter δ.

.
3.2. Numerical simulations for condensation transition at steady state

In order to verify the theoretical predictions of the scaling laws by previous GCE analysis, we perform the following Monte Carlo (MC) simulations. The underlying SF network is generated using the Barábasi–Albert model,[23] with a degree distribution with γ = 3, and hence leads to a critical value δc = 0.5. The initial conditions for the simulations are random distribution. To ensure that the sampled data is taken from the steady state, we keep the system running extra time after the realization of the steady state, and then measure the quantities interested.

To see the presence of the finite density effect, in Fig. 1, the scaling relations of mk versus kρδ are shown at δ = 0.2 < δc for different network sizes L = 1000, 2000, 4000, and 8000, and various densities, ρ = 0.5 (purple open diamonds), ρ = 2 (blue open triangles), ρ = 5 (red open squares), ρ = 10 (green open circles), respectively. As expected from Eq. (33), for k < kc (i.e., on the left hand side of the crossover degree kc), it is obvious that mk ∼ kρδ Lδ – δc, hence for given network with size L, mk ∝ kρδ for a fixed δ; while for k > kc (i.e. the right hand side of the crossover degree kc), mk ∼ k1/δ ρ L1 – δc/δ, for a given network with size L, mk ∼ k1/δ ρ ∝ (k ρδ)1/δ. In other words, in both the cases, the average occupation number is an scaling function of kρδ at given L, namely, , where is some scaling function. One finds that in Fig. 1 the simulation data nicely collapse onto the same curve, hence fully support our theoretical predictions.

Fig. 1. The scaling relations of the average occupation number mk versus kρδ at δ = 0.2 for different SF network sizes: (a) L = 1000, (b) L = 2000, (c) L = 4000, (d) L = 8000, at various densities ρ = 0.5 (purple open diamonds), ρ = 2 (blue open triangles), ρ = 10 (red open squares), ρ = 50 (green open circles), respectively.

The predicted finite density effects can be viewed from the other perspective. From Eqs. (23) and (31), a ρ-dependent ratio can be derived, i.e., kc(ρ2)/kc(ρ1) = (ρ1/ρ2)δ. This ratio eliminates the effect of the proportional constants in kc. We must emphasize at this point that, without taking into consideration of the variable density, this ratio should be a density independent constant, i.e., kc(ρ2)/kc(ρ1) = 1. In Fig. 2 we use kc at ρ1 = 2 as a benchmark, plotted theoretical ratios kc(ρ)/kc(ρ = 2) with various densities, ρ = 1 (black solid line), ρ = 5 (blue solid line), and ρ = 10 (red solid line), respectively, as a function of δ at given network size L = 4000. The purple dashed horizontal line shown in Fig. 2 just corresponds to the density independent ratio, which is unity obviously. However, MC simulations show the density dependent ratios for ρ = 1 (open black diamonds), ρ = 5 (open blue triangles), and ρ = 10 (open red squares), respectively, and agree with our theoretical predictions very well. In principle, kc’s in δ > δc regime do not exist due to the existence of lower bound kmin for densities ρ = 5 and ρ = 10, hence we formally draw the theoretical results by blue (ρ = 5) and red (ρ = 10) dashed lines.

Fig. 2. The presence of the finite density effect in the crossover degree ratios kc(ρ)/kc(ρ = 2) versus the parameter δ with an SF network size L = 4000, for various densities ρ = 1 (black open diamonds and solid line), 5 (blue open triangles and solid line), and 10 (red open squares and solid line), respectively. The solid lines correspond to the theoretical calculation and the symbols to the MC simulation data. Note that kc’s do not exist in regime δ > δc (δc = 0.5 is the vertical red dotted dash line) for ρ = 5 and ρ = 10 due to the lower bound of degree kmin, two theoretical dashed lines are only formally drawn in that region.
4. Finite density scaling law of relaxation dynamics

The preceding theoretical predictions show that the occupation number distribution in the steady state is determined not only by the degree distribution being independent of any other characteristics of underlying networks, but also the density of the particles hopping on the networks. However, the relaxation dynamics of ZRP is getting more complicated due to the effects of both the structure of underlying networks and the density.

So far only has a little been known about the dynamics of condensation in ZRP.[17,37,39] The finite density effects on the relaxation dynamics is hence one of the important issues to be addressed. In the following we only consider the case of δ < δc, in which a condensation phase transition occurs. We rely on the MC simulations and the scaling ansatz in order to understand the dynamical scaling properties.

Scaling ansatz is helpful for obtaining the qualitative features of relaxation dynamics. Previous studies[17,37] revealed that particles form macroscopic condensate by successive coarsening processes of the small clusters. The smaller clusters merge into larger ones, and grow until a single macroscopic condensate forms. The scaling hypothesis leads to the power-law growth of the relaxation time with the system size as tR ∼ Lβ, where β = 1 – δ is the dynamic exponent.

The highly heterogeneous structure of SF networks deeply affect the relaxation dynamics of interacting particles system defined on them. It was hypothesized in Ref. [27] that the relaxation process follows some hierarchical dynamics. There exist two characteristic time-dependent degree scales kΩc(t) and kΩ(t), which play the role of the crossover degree kc and the role of kmax, respectively. During transient time t, a subnetwork of small degree nodes is equilibrated first and the equilibrated subnetwork keep expanding until kΩ(t) grows in time and eventually reaches kmax.

It is natural to assume that a similar hierarchical dynamics occurs in a variable finite density situation. Starting from a fully random distribution at initial time, all particles behave as random walkers without interaction in a short period of time. Later on the distribution of particles further evolves and eventually approaches the steady state. Suppose that there exist two characteristic degree scales in hierarchical relaxation dynamics. Those nodes with degree k ≤ kΩ consist of a smaller equilibrated subnetwork, denoted by, say, Ω, of size LΩ < L in a transient time t ≪ tR, where tR is the relaxation time. Inside this subnetwork, the other characteristic degree kΩc plays the role of crossover degree kc of the whole network. The equilibrated subnetwork keeps growing to proceed to the higher hierarchy until kΩ reaches kmax of the network, the steady state of the ZRP is then formed. We may refer to those nodes with kΩ as the temporary hub.[28] The average occupation number mkΩ versus the transient time t, according to the hierarchical dynamics, scales as

with ν = 1 – δ being the dynamic exponent.[28] Clearly kΩ evolves with time, however, we emphasize at this moment that kΩ is also density dependent, i.e.,

We now turn to the other scale in the transient state: the crossover degree kΩc. Just like kc’s scaling in the steady state of the network, it scales as

and again, a density dependent factor ρ−δ comes in.

The average occupation number of a node within the equilibrated network mk reads

As a result, by applying Eq. (35), the scaling law of time- and density-dependent mk can be obtained as

These explicit scaling laws of mk(t) versus kρδc for some transient time t can be unified in a single relation , in which is some scaling function, for given t and δ. MC simulations have to be performed to verify the above analytical predictions.

The results of MC simulations on the finite density effects of the relaxation dynamics are shown in Fig. 3(a). Note that the different time scales are involved in our relaxation dynamics simulations: the transient time scale t, the relaxation time scale tR, the evolution time scale of the system τ, and t ≪ tR ≪ τ. We took the evolution time τ = 220, which is long enough for system to evolve towards the steady state. The average occupation number mk(t) at transient time t = 210, versus kρδc with a given network size L = 4000, but with various densities ρ = 0.5 (purple open diamonds), ρ = 2 (blue open triangles), ρ = 10 (red open squares), and ρ = 50 (green open circles), respectively, are shown.

Fig. 3. Data collapse plot of (a) the average occupation number mk versus kρδc at the transient time t = 210, (b) the evolution of ratio of IPR to its steady state limit, It/I∞, versus t/ρν, with ν = 1 – δ being the dynamic exponent in the relaxation dynamics, at δ = 0.2 for network size L = 4000, and various densities ρ = 0.5 (purple open diamonds), ρ = 2 (blue open triangles), ρ = 10 (red open squares), and ρ = 50 (green open circles), respectively.

As one can see from the left panel in the figure, the collapse of simulation data of mk w.r.t. kρδc indicates the scaling law in the evolution of the relaxation for different densities. To understand the condensation dynamics, an appropriate characteristic order parameter is the infinity time limit of inverse participation ratio (IPR) It,[27] which is defined as

According to the hierarchical dynamics ansatz in condensation processes, the occupation number of particles on the temporary hub dominates in an equilibrated sub-network, note that the degree of the temporary hub kΩ scales as tδc/ν/ρδc as in Eq. (35), combining with , the size of the equilibrated sub-network satisfies LΩ ∼ t/ρν for a given transient time t. Hence the IPR ratio, describing the extent of condensation for a given network size L, possesses the same scaling relation It/I∞ ∼ t/ρν, while for the long run τ ≫ tR the steady state is reached and the particles on the hub dominate.

In Fig. 3(b), we plot the scaling relations of It/I∞ versus t/ρν at δ = 0.2 for ρ = 0.5, ρ = 2, ρ = 10, and ρ = 50, respectively. A total τ = 220 time-step evolution has been fulfilled to ensure the realization of the steady state, i.e., Iτ → I∞. One immediately sees that the increasing density naturally postpones the occurrence of the condensation transition, and the collapse of data of different densities nicely confirms our prediction.

5. Conclusions

In conclusion, we have explored the finite density effects on both the steady and dynamical properties of the ZRP on SF network with hopping rate function u(n) = nδ. By means of the GCE method, we find that the steady state condensation phase transition is driven not only by the disorder of the underlying SF network, but also by the finite density of the interacting particles. The crossover degree kc and the average occupation number mk are found to be both density dependent, and satisfy the corresponding scaling laws. In contrast to the ZRP on SF network with fixed density, the ratios of the crossover degree kc for two different densities at given size of network exhibit a non-constant behavior.

The influences of the density on the relaxation dynamics is also analytically investigated. The hierarchical characteristics of the relaxation dynamics renders the process proceeding from nodes with lower degrees toward those with higher degrees. With the help of the scaling ansatz, the scaling laws of the condensation transition (to the steady state) are predicted. At a specific time t ≪ tR the average occupation number scales as mk ∼ kρδc for various densities. The evolution of the ratio of IPR to its steady state limit follows a scaling law It/I∞ ∼ t/ρν with dynamical exponent ν = 1 – δ. To verify our analytically derived scaling laws, we have performed the Monte Carlo simulations for varied densities. The corresponding scaling relations are validated very well by numerical experiments.

Reference
[1] Anderson M H Ensher J R Matthews M R Wieman C E Cornell E A 1995 Science 269 198
[2] Chowdhury D Santen L Schadschneider A 2000 Phys. Rep. 329 199
[3] van der Meer D van der Weele K Lohse D 2002 Phys. Rev. Lett. 88 174302
[4] van der Meer D van der Weele K Lohse D 2004 J. Stat. Mech.: Theor. Exp. 04 P04004
[5] Krapivsky P L Redner S Leyvraz F 2000 Phys. Rev. Lett. 85 4629
[6] Bianconi G Barabási A L 2001 Phys. Rev. Lett. 86 5632
[7] Su G F Zhang X B Zhang Y 2012 Eurphys. Lett. 100 38003
[8] Pathria R K Beale P D 2011 Statistical Mechanics 3 New York Academic Press
[9] Liggett T M 1999 Stochastic Interacting Systems: Contact, Voter, and Exclusion Processes Berlin Springer
[10] Schmittmann B Zia R K P 1995 Statistical Mechanics of Driven Diffusive Systems Domb C Lebowitz J L New York Academic Press
[11] Evans M R Kafri Y Koduvely H M Mukamel D 1998 Phys. Rev. Lett. 80 425
[12] Evans M R Hanney T Majumdar S N 2006 Phys. Rev. Lett. 97 010602
[13] Spitzer F 1970 Adv. Math. 5 246
[14] Evans M R Hanney T 2005 J. Phys. A: Math. Gen. 38 R195
[15] Godrèche C 2007 Lect. Notes Phys. 716 261
[16] Evans M R 2000 Braz. J. Phys. 30 42
[17] Großkinsky S Schütz G M Spohn H 2003 J. Stat. Mech.: Theor. Exp. 113 389
[18] Majumdar S N Evans M R Zia R K P 2005 Phys. Rev. Lett. 94 180601
[19] Godrèche C Luck J M 2012 J. Stat. Mech.: Theor. Exp. 2012 P12013
[20] Albert R Barábasi A L 2002 Rev. Mod. Phys. 74 47
[21] Eguiluz V M Chialvo D R Cecchi G A Baliki M Apkarian A V 2005 Phys. Rev. Lett. 94 018102
[22] Boccaletti S Latora V Moreno Y Chavez M Hwang D U 2006 Phys. Rep. 424 175
[23] Barabási A L Albert R 1999 Science 286 509
[24] Dorogovtsev S N Goltsev A V Mendes J F F 2008 Rev. Mod. Phys. 80 1275
[25] Barthelemy M Barrat A Pastor-Satorras R Vespignani A 2004 Phys. Rev. Lett. 92 178701
[26] Sood V Redner S 2005 Phys. Rev. Lett. 94 178701
[27] Noh J D Shim G M Lee H 2005 Phys. Rev. Lett. 94 198701
[28] Noh J D 2005 Phys. Rev. 72 056123
[29] Tang M Liu Z Zhou J 2006 Phys. Rev. 74 036101
[30] Cohen R Erez K ben-Avraham D Havlin S 2000 Phys. Rev. Lett. 85 4626
[31] Noh J D Rieger H 2004 Phys. Rev. Lett. 92 118701
[32] Krug J 2000 Braz. J. Phys. 30 97
[33] Evans M R 1996 Europhys. Lett. 36 13
[34] Jain K Barma M 2003 Phys. Rev. Lett. 91 135701
[35] Evans M R Majumdar S N Zia R K P 2004 J. Phys. 37 L275
[36] Zia R K P Evans M R Majumdar S N 2004 J. Stat. Mech.: Theor. Exp. 2004 L10001
[37] Godrèche C 2003 J. Phys. A: Math. Theor. 36 6313
[38] Großkinsky S Hanney T 2005 Phys. Rev. 72 016129
[39] Godrèche C Drouffe J M 2017 J. Phys. A: Math. Theor. 50 015005